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1 Summary 


In this supplemental PDF we: 

e Provide expanded plots and description of our comparison of ST rough dielectric BSDFs to measured 
data (subsection 2.1) 

e Give a detailed derivation of the Beckmann-superposition vNDF sampling scheme for the Student-T 
NDF (section 2) 

e Discuss our null-collision algorithm for sampling and evaluating multiple scattering from rough 
and porous generalized Smith surfaces with NDFs extending to the full sphere, together with some 
new parametric NDFs and their properties (section 3) 

e Discuss details of our BSDF interface that permits easily implementation of new null and shape 
invariant NDFs (section 4) 

e Discuss several results that were created using the above tools (section 5) 


2 Student-T NDF: New Sampling and Motivation 


2.1 ST comparisons to Rough Glass 


The GGX NDF (equivalent to the Trowbridge-Reitz NDF) was made popular in graphics by Walter et al. 
(2007) after using it to approximate the behaviour of ground glass, noting that Beckmann’s NDF was a 
less accurate choice. In the same work, it was observed that several glass samples (etched and frosted) 
appeared to exhibit a behaviour somewhat between GGX and Beckmann, with neither constituting an 
optimal choice. 

The Student-T (ST) NDF (Ribardière et al., 2017) was later introduced and conveniently interpolates 
between Beckmann (y = œ) and GGX (y = 2) via its dependence on an additional shape parameter y > 3/2. 
The primary reasons for proposing this NDF were convenient mathematical properties that are required 
for its efficient use (such as closed-form shadowing/masking factors). To the best of our knowledge, ST 
BSDFs have not been compared to measured data in order to determine if ST is indeed a good physically- 
motivated choice for interpolating between Beckmann and GGX. 

In this section, we revisit the original data from Walter et al. (2007) and compare predictions of the ST 
microfacet BSDFs to the measurements for the ground, etched, and frosted glass samples. In Figures 1 - 3, 
we plot each measurement from Walter et al. (2007) and compare to Beckmann, GGX and two ST BSDFs: 
one where y = 3 is fixed and one where y is free to vary for that sample. Our plots differ from those 
of Walter et al. (2007) by plotting VL to provide an approximate gamma correction, which allows closer 
inspection of the errors in the darker regions of the measurements. In the case of frosted and etched glass, 
we see that the ST NDF not only outperforms GGX and Beckmann, but in fact matches the measurements 
of both samples as accurately (or more) than GGX matches the ground glass measurements. From this we 
conclude that ST is a very strong choice for a two-parameter model that can blend between the statistics 
of Beckmann and GGX surfaces. Renders of rough glass using these settings and comparisons to Walter’s 
original GGX and Beckmann fittings are provided in the supplemental images folder. 

To put the ST NDF on equal computational footing to Beckmann and GGxX, in the next section we 
derive vNDF sampling, which is required for low-variance single-scattering sampling, as well as reference- 
quality random walk multiple-scattering (required to avoid large energy loss at high roughnesses as seen 
in (Ribardière et al., 2017, Fig.7)). 
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(a) slightly tweaked 2007 fit with Beckmann 


Frosted - ST3 fit - Single-scatter only, scale = 0.861, a=0.454 


(c) Student-T with y = 3 


Frosted - GGX fit - Single-scatter only, scale = 0.861, a=0.454 


(b) 2007 fit with GGX 


Frosted - ST fit - Single-scatter only, scale = 0.752417, a=0.428006, y=2.80965 


(d) general Student-T fit 


Figure 1: Frosted glass: note how for the two most grazing measurements in each plot (peaking around 
220 and 240 degrees) the ST fits (c-d) are more accurate than previous models (a-b). Note also the improved 
accuracy of ST for the normal (180) profile for low intensities. 


Etched - Beckmann fit - Single-scatter only, scale = 0.711, a=0.493 
20F r r T T r r 7 J 


(a) 2007 fit with Beckmann 


Etched - ST3 fit - Single-scatter only, scale = 0.835372, a=0.52379 


(c) Student-T with y = 3 


Etched - GGx fit - Single-scatter only, scale = 0.955, a=0.553 


2.0F T T T T r r T 7 
. 
41.5} J 
Bart 
. 
1.0} ] 
cme 
0.5 ae 
i fe 
Se, 
b. "ee 
*eecoss cep at 
0.0 HLI Í oe ooo sses sce ae 
: 1 
100 120 140 160 180 200 220 240 
8o 


(b) 2007 fit with GGX 


Etched - ST fit - Single-scatter only, scale = 0.840502, a=0.523738, y=2.82261 
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(d) general Student-T fit 


Figure 2: Etched glass: note a similar improvement of ST over Beckmann or GGX for both low intensities 


and for the more grazing (220, 240 peaked) profiles. 


Ground - Beckmann fit - Single-scatter only, scale = 0.542, a=0.344 Ground - GGxX fit - Single-scatter only, scale = 0.755, a=0.394 
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(a) 2007 fit with Beckmann (b) 2007 fit with GGX 
Ground - ST3 fit - Single-scatter only, scale = 0.678133, a=0.381009 Ground - ST fit - Single-scatter only, scale = 0.805458, a=0.405877, y=1.79852 
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(c) Student-T with y = 3 (d) general Student-T fit 


Figure 3: Ground glass: in this case the GGX distribution (b), which was motivated by this exact physical 
sample, performs best and is not improved upon by Student-T. 


2.2 Review 


Before deriving our ST sampling procedure, we begin with some necessary review. With notation u = 
cos 0, the Beckmann NDF with roughness a is (Blinn, 1977) 


D®(u, a) = —~0lu), (1) 


where ©@(u) is Heaviside’s function. The GGX NDF is 
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and the student-T NDF is 
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The corresponding slope distribution is 
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with marginal slope distribution 
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The scattering cross-section is 
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2.3 Student-T as a superposition of Beckmann surfaces 


The key insight to our method is to note that the ST NDF can be expressed as a gamma-superposition of 


Beckmann NDFs: 
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We will use the method of continuous decomposition ((Butler, 1954, p.254), (d’Eon, 2022, Sec. 70.3.3)) 
for sampling non-uniform random variates, where a target distribution is expressed as a superposition of 
another parametric distribution, for which a sampling procedure is known. Here the parameter m > 0 
indexes the superposition, and for each m the roughness of the corresponding Beckmann distribution in 


the mixture is æy (y — 1)/m. 


2.4 Sampling the visible distribution of normals 


For a ray arriving straight down along the surface normal, we could sample the ST vNDF in a two stage 
procedure using Equation 7 to first sample m from the gamma distribution 
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(using Marsaglia’s method or C++ 11 code) and then sample the Beckmann NDF with roughness ag = 
av (y — 1)/m. However, when the incoming direction w; is not along the normal, we require a different 
distribution, due to how the various Beckmann distributions in the continuous mixture self shadow each 
other. We can find this distribution by repeating the slope-space derivation of Heitz and d’Eon (2014), but 
with Equation 7 as the NDF. We find that we need to sample slope p from the marginal slope distribution 
(Equation 5) 


(8) 


u/V1-u2 
/ P3"(p,a,y)(u — pV1 — u?)dp, (9) 


o0 


which we expand using the gamma superposition (Equation 7) to find 
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Combining the last two equations, and integrating over p, we are left with a distribution over m. Using 


stretch invariance, we can solve the problem for roughness @ = 1 and linearly scale the sampled Beckmann 
roughness by a. We find for a = 1, 
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gives the probability density for the Beckmann NDF indexed by m to contribute to the ST vNDF for a ray 
arriving along cosine u. Its normalization is simply the cross section for the ST NDF (Equation 6), 


f fo u)dm = o(u) = (1+ A(u))u. (13) 
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If we can sample m from F (m, u) given u, then we can then sample the Beckmann vNDF with roughness 
ap = ay (y — 1)/m and we have sampled the ST vNDF. 
2.5 Sampling m from the upper hemisphere 


For u > 0, the density f in Equation 11 is a sum of three non-negative densities, each involving gamma 
distributions. We can sample the sum by discrete decomposition (Everett and Cashwell, 1983, C3. p.55). 
The term involving erf can be sampled using a rejection procedure. We first decompose f into three terms 
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where f; are normalized densities over m € [0, co]. We define 
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We find (by integrating each over m > 0 in Mathematica) that the discrete selection probabilities for the 
three densities are 
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The final procedure for sampling the visible distribution of normals for the ST NDF for an incident direc- 
tion pointing down towards the surface with cosine 0 < u < 1is to first select a lobe f; using discrete prob- 
abilities {pi(u), p2(u), p3(u)}, then sample m from the corresponding density {fı (m, u), fo(m, u), fa(m, u)}, 
and finally sample the Beckmann vNDF with roughness ag = æy (y — 1)/m. 

Densities fı and fz are gamma distributions and can be directly sampled using C++ 11 methods. Specif- 
ically, fı is a gamma distribution with shape parameter a = y — 3 and scale parameter b defined above. fz 
is a gamma distribution with shape parameter a = y — 1 and shape parameter b = 1. We finally note that fz 
is a gamma distribution with shape parameter a = y — 1 and scale b = 1 that is modulated by the function 


erf (} =} which can be sampling using rejection by repeatedly sampling m from the corre- 


sponding gamma distribution and rejecting and trying again each time &; < erf g e) where 


& € [0,1] are uniform variates. We observed the expected number of gamma random variates required 
to sample the vNDF is 1.2, so the rejection of f; happens infrequently. We remark that for y € {2, 3, 4} it 
is possible to exactly sample fz using a fixed number of random variates. 

The probabilities p, and p, are inconvenient functions of cosine u and include the hypergeometric 
function 2F1. For these reasons, we derived (using TuringBot) approximate equations (provided in the 
supplemental code) for each, which we observed to work very well in practice. 


2.6 Upwelling case 


The down-welling case is sufficient for BSDF sampling the single-scattering portion of the rough surface. 
However, for a full multiple-scatter random walk, we also need to handle vNDF sampling for up-welling 


rays (—1 < u < 0) (Heitz et al., 2016). Unfortunately, this requires a different sampling procedure because 
the three terms in the m-sampling density f(m, u) are no longer non-negative, excluding the possibility 
of sampling by discrete decomposition. 

We propose an approximate solution using generalized gamma distributions. We observed that the 
target distribution f(m,u < 0) is well approximated by a generalized gamma distribution with three 
parameters a, b, and c, with probability density 
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This distribution can be sampled by sampling a gamma variate x with parameter a and width b = 1, 
and then returning m = bx'/° (Tadikamalla, 1979). We used a fitting procedure (TuringBot) to find the 
constants a, b, and c, given cosine u and shape y. We determine b first by requiring that it match the exact 
first moment of the target distribution, 
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Our fitting experiments found approximations for a and c, 
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2.7 Discussion 


We have tested the above sampling procedure and found excellent agreement compared to a sampling of 
the ST vNDF using a null collision (rejection) procedure and verified that the supplemental code produces 
close results to GGX and Beckmann BSDF values when y is appropriately set. 

This sampling scheme has several limitations, such as its use of a theoretically unbounded number of 
random variates per sample, as well as the requirement for gamma variates. In practice, we observed the 
rejection scheme required for the down-welling case to happen about 20% of the time. A further limita- 
tion is that both the superposition sampling and rejection will hinder the application of low discrepancy 
sampling schemes. 


3 A Null-Scattering Formulation of Rough Surface Scattering 


In this section we derive a new Monte Carlo procedure for stochastically evaluating and sampling Smith 
microfacet BRDFs. This approach is motivated by two goals: 1) to allow quick and easy addition of new 
NDFs to a material system, and 2) to expand on the proposal of Dupuy et al. (2016) to formulate a Smith 
theory of microfacet scattering for highly rough materials where some microfacets face into the downward 
hemisphere, thereby violating the height-field assumption of the standard microfacet theory (Smith, 1967). 
Smith’s theory is a powerful general method for modeling surface scattering and supports anisotropy 
(Heitz, 2014), importance sampling (Heitz and d’Eon, 2014) and multiple-scattering (Heitz et al., 2016) for 
many surface statistics (Walter et al., 2007; Ribardière et al., 2017). However, the height-field assumption is 
not universally appropriate, and deriving Lambda functions and vNDF sampling procedures is not always 
possible. 

Dupuy et al. (2016) describe a generalization of the Smith theory for a volume where the NDF extends to 
support the full sphere of directions. This was accomplished by noting that all Smith heightfield volumetric 
models correspond exactly to a generalization of a uniform anisotropic/microflake (Jakob et al., 2010) half 
space of spatially-independent scattering particles/flakes. The total anisotropic macroscopic cross section 
0;(w) governing this volume was found to be the normalization factor of the distribution of visible normals 
for the given NDF D(w) of the material, via the relationship (Dupuy et al., 2016, Eq.(6)) 
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where (@ ,@2) > 0 is a clamped dot product. The volumetric equivalence in Equation 23 immediately 
and uniquely determines the required random walk: given any NDF D(wm) > 0 the free-path length 
distribution for sampling random walks in the volume along direction w is then 


p(s; 0) = (wes, (24) 


After each collision (moving in direction w), a microfacet normal wm is then sampled from the normal- 
ized distribution Dyjs(@m, @) and the BSDF of that microfacet is sampled to determine the new direction. 
This process is continued until the ray escapes the half space. However, this requires accurate numerical 
evaluation of o, and importance sampling of Dyis(@m, ©). 

No specific NDF or related procedure for sampling distances and visible normals in such a medium 
was presented by Dupuy et al. (2016). In a followup paper (d’Eon, 2016) it was shown that the cross section 
for the von-Mises Fisher (vMF) distribution was complicated, and no sampling procedure for the related 
distribution Dyis(@m, @) has been found. To the best of our knowledge, no practical numerical scheme for 
evaluating a non-uniform full-sphere NDF BRDF has been presented. 

In this talk, we show how an approach inspired from Woodcock/delta tracking for sampling distances 
in inhomogeneous volumes (Novak et al., 2018) can stochastically evaluate any generalized Smith model 
given only the NDF and its majorant as input, requiring no additional derivations. Where delta-tracking 
in inhomogeneous volumes spatially homogenizes the cross section o;(x), we use an analagous proce- 
dure to homogenize the cross section in the angular domain o;(w) (it is already spatially homogeneous). 
This is achieved by introducing just enough virtual microflakes to homogenize the cross section to be a 
constant over the sphere of microfacet normal directions. Collisions with the new medium result in some 
(directionally-dependent) portion of those collisions to be virtual/null, resulting in the rejection of these 
collisions, with the random walk continuing in the same direction as if the null flake was never present. 


3.1 The Algorithm 


We consider a homogeneous half space of scattering and absorbing particles/flakes. Spatial homogeneity 
implies no dependence on position for the cross section o;(w), but the angular dependence via direction w 
is required to exhibit the correlations of surface scattering, with (in general) asymmetry o;(@) # o;(—@) 
(Dupuy et al., 2016). To sample the microfacet BRDF, given incoming direction w; and general NDF, 
possibly defined over the full sphere, we use the following algorithm: 

e Given an NDF D(@m) > 0, we define the constant majorant NDF D(@m) > D(@m) whose value 
is the maximum that the NDF D takes on (typically the value straight up along the macrosurface 
normal). 

e We stochastically evaluate the microfacet BRDF by sampling a random walk in the half space, using 
the homogeneized volumetric Smith model Dupuy et al. (2016). The walk begins at the boundary 
along with initial direction œw = —w;, where «; is the common BRDF notation of a vector pointing 
away from the surface. 

e We sample all collisions in the medium using null scattering. Specifically, we collide with the an- 
gularly homogeneized distribution of microfacets (with NDF D(@m)). Combining Equation 23 with 
the uniformity of the NDF results in o; is being a constant. Without loss of generality, since we are 
computing a BRDF with no lateral displacement (sometimes known as “position-free simulation” 
Guo et al. (2018), even though the depth-tracking is necessary), and because the surface is expanded 
into a half space, we can assume o; = 1. 

e While no actual particles are traced against, the following geometrical procedure illustrates the 
core sampling idea. For each sampled collision (real or ficticious), we consider a small spherical 
particle with the pair of NDF values (real and ficticious) mapped to its surface (Figure 4). We then 
consider a uniform field of rays along direction w intersecting this sphere, and sample from this set 
of rays uniformly at random (using disk sampling following by rotation to w and then a ray-sphere 
intersection test). This produces a microfacet normal wm that is the surface normal of the sphere 
at the intersected location. Using this microfacet normal direction we consider that this sampled 
collision is real with probability D(@)/D(@m) and otherwise ignore the collision and sample a new 
distance with w unchanged, analogous to delta tracking in inhomogeneous media (see Figure 4). 

e For each real collision, the BSDF of the microfacets is sampled using w and w,, to determine the 
next direction of the random walk. 

e Upon exiting the half space, the particle weight (which may deviate from 1 for absorbing micro- 
facets) and direction w are returned. 

This procedure is equivalent to Smith’s model in the case that the NDF is restricted to the upper 

hemisphere. We have verified this via MC simulation for Beckmann and GGX NDFs. 

An important property of this algorithm is that it jointly samples free-path-lengths p(s; œw) and vis- 

ible normals Dyis(@m) in a single step, without having to know o;(w). The disk projection of uniform 
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Figure 4: Illustration of our stochastic null-scattering formulation of non-heightfield Smith scattering. We 
model a rough surface as a homogeneous half space of scattering particles with no spatial correlation. 
Each time a particle collision is sampled, we select a random ray along the current direction (thick line) 
from a uniform field of rays that intersect a virtual spherical particle (thin lines). At the intersection point 
of the sampled ray with the particle (black dots), we consider the amounts of real and ficticious density in 
the NDF at the sampled microfacet direction. Here the real and ficticious NDFs are plotted parametrically, 
offset from the particle to illustrate their relative magnitudes. In the illustrated example here (the sampled 
thick path), roughly 60% of the total density of the uniform majorant NDF in the sampled direction is 
real, and so the particle is stochastically determined to cause a true collision with probability 60%, and 
the sampled microfacet normal is used with the BRDF of the microfacet to sample the next direction. This 
process continues until the path exits the half space, and the final direction and weight of the path are 
returned. 


rays onto the sphere exactly accounts for the geometry term (the clamped dot product) in Equation 23 
and the rejection procedure exactly normalizes the true distribution, and by Equation 23 produces the 
desired exponential mean free path length o;(w) (this follows from properties of Poisson processes, see 
e.g. Georgiev et al. (2019)). 

Another key feature of this algorithm is that only D(w,,) is required to implement sample(). This 
permits adding almost any (bounded) NDF to a renderer. 


3.2 Next-Event Estimation 


In order to implement eval(), we require the same random walk as in the sample() method, but with NEE 
at each real collision along outgoing direction w,. This requires knowing the extinction in direction wo, 
which depends on the cross section o(wo) of the NDF. In practice, this function is unavailable in closed 
form for many NDFs. To circumvent this limitation, and simplify the addition of new NDFs to the material 
library, we use the cross section of the Dirac NDF (described next) as a Green’s function and numerically 
compute o(@,) using a 1D quadrature over the cosine of the microfacet orientations (note that this step 
of the algorithm is limited to symmetric NDFs that depend only on direction cosine u). 


3.2.1 The Dirac NDF 


The general case of a Dirac ring on the sphere of normal directions with azimuthal symmetry corresponds 
to a surface where all normals inclined by a common cosine, uy. This might arise in the case of a surface 
that is thoroughly pitted by the tip of a hard cone oriented along the mean normal (Figure 5). The Dirac 
NDF is simply 

D(u) = (u — un), (25) 


and leads to a piecewise simple expression for the cross-section (d’Eon, 2016) 
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Figure 5: A smooth gold surface pitted with a conical tip aligned to the surface normal creates a micro- 
surface with an isotropic NDF given by a Dirac delta of the cosine of the microfacets to the macrosurface 
normal. 


Family of Dirac cross sections for azimuthally-symmetric normal distributions 
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Figure 6: The Dirac-delta cross sections os(u, un) for a variety of u, form a family of normalized functions 
blending between z(u + |u|) and 2V1 — u?. 


The Dirac NDF cross section can be decomposed into three functions, 
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We can also compactly write 
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functions (Figure 6). The Dirac cross section is closely related to a cross section for random plates in 
anisotropic random media studied in tree canopies (Nilson, 1968b,a; Myneni et al., 1988; Knyazikhin and 
Marshak, 1991). Similar expressions to Equation 26 appear in (Davison, 1957, p.235,Eq.(17.12)). 

We can use the Dirac cross section as a Green’s function and express the cross section of any symmetric 
normal distribution D(u) as, 


o(u) = ie o3(u, u’)D(u')du’. (33) 
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To acquire o(u) during rendering having only implemented D(u), we numerically integrate Equation 33 
using a Gaussian quadrature. To be conservative and non-adaptive, we chose 100 points, but fewer is 
likely sufficient in many cases, especially given the smooth NDF shapes likely to be efficient with the 
null-scattering algorithm. 


3.3 Limitations 


While the null-scattering approach for general microfacet BSDFs is very flexible, it has some important 
limitations. Highly-peaked (low roughness) NDFs will have very high rejection rates for almost all di- 
rections, not just upwelling ones, causing the method to perform poorly. In our experience, roughnesses 
below o = 0.4 are not recommended. 

Path tracing typically also requires a deterministic pdf() (if only approximate), in order to use MIS in 
a unidirectional path tracer, for example. For our tests using this method in Mitsuba 0.6 (Jakob, 2010), we 
used a Lambertian pdf in all cases. Given that the approach is best suited to highly rough materials (which 
have NDFs on the full sphere), we found that this compromise was adequate for an initial exploration of 
the space of biscale and porous BSDFs. 

Both of these limitations leave open important areas for future work. 


3.4 More NDFs 


With the ease of adding new NDFs to the library, we note several additional parametric families that do 
not seem to have gained attention previously in rendering and might find use in some applications. 


3.4.1 Bessel K NDF 


We define a parametric NDF that interpolates the exponential (defined below) and Beckmann NDFs. 
Specifically, we define the Bessel K shape invariant NDF with roughness a and shape parameter a > 0 
using continuous superpositions of Beckmanns to find 


œ s-a’ yy-2 
DS = C A DE (u, ava) da’ 34 
(u, a, a) f PGi (ua a’) a (34) 
co g1-ag-% g’?a=1 z 
= ny 2) "1N2) da’ 3 
i Ta (u. aa 2) a (35) 
2(1- u2) 2-2 (au) K-a (ac in) 
~ maul (a) ee 


where K,(x) is the Bessel-K function. The corresponding slope distribution is 


PX (p,q, æ, a) = 37 
22 (P, q, &, a) nT (a) (37) 
with marginal slope distribution 
2a-*- 3 |pl@-2K,_, (22) 
PS ,&,a) = 38 
2 (p, æ, a) TTO (38) 
The Smith A function for the Bessel K NDF can be written 
o u piai ia a (ga) 
A` (u, a, a) = = + —— 
2 VaT (a) 
~ 2 
L ( ) yru iF, ER 3 -a - Gz) ar a 2-4 F ( 1 i u2 
-yru sec( zra a tu —u 1Fo |a;a+ —, a + 1; -— 
2 avi -uT (a) 2 (u? — 1) @? 
(39) 


For integer values of a, the A functions reduce to simpler expressions. For half-integer values of a, Meijer 
G functions appear (see BesselK.nb for more details). 

As a — œ, the NDF DX (u, @/-Va,a) approaches the Beckmann NDF with roughness a (Bahar and 
Fitzwater, 1983). For a = 3/2, the Bessel K NDF reduces to the exponential NDF (below). Our definition is 
a generalization of a slope distribution proposed by Bahar and Fitzwater (1983) for integer values of a (to 
follow their constant-mean squared slope parameterization use DÝ (u, æ / ya, a)). 
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The Bessel Ko NDF (a= 1) For this value of a the resulting NDF is singular at u = 1 and has explicit A 
and cross section: 


2Ko (are ae) 


DŽ (u, æ, 0) = ——O(1 - u) (40) 
mou 
_ 2u 
V1—- 2 aVi-u2 
AK(u,a,0) = Č "e , 0<u<1 (41) 
4u 
_ 2u 
1 2u V1 — ue «Vizu? 
oX(u,a,0) = 701 — u? O(—u)e«V1-u2 + uO(u) 2 = +1 (42) 
u 


The Ko NDF was also considered by Bourlier et al. (2002). 


The Exponential NDF (a = 3/2) Microsufaces with an exponential distribution of slopes were pro- 
posed by several authors (Hughes, 1962; Muhleman, 1964; Hagfors, 1966; Barrick, 1968). Using a different 
parameterization of roughness, we define the the exponential shape invariant NDF to be 


© 96-0 far _ 2Ni-u2 
DE? (u, æ) = £ NT DB(u, aVa’)da’ = æ = ei - u). (43) 
0 va Tau 
The corresponding slope distribution is 
_ 2V pg 
pË a 
Pi" (p,q) = s (44) 
with marginal slope distribution 
2lpl 
ies alpiK; (22) 
P; (p, a) = ~ (45) 
The Smith A function for the Exponential NDF can be written 
1 
oi etal 40) 2 (ai) 
1-w)at! 1 3 2uK: 
APP (u, a) = iii ak aL) A Al : =f (46) 
T mTavi—u2 2 


in terms of Meijer G functions, which has an equivalent representation using Struve L functions Brown 
(1980). Bourlier et al. (2002) also considered Smith shadowing functions for exponential distributions of 
slopes. For accurate approximations of the Exponential NDF Lambda function, see Mathematica source 
code at: ExponentialNDF.nb. 


4 Facet Forge 


In the supplemental codebase facet-forge, we provide a C++ implementation of the above BSDF tools 
using a new object-oriented interface that minimizes the amount of code required to add a new NDF. By 
deriving from a common base NDF class, both a ShapeInvariantNDF and a NullNDF are implemented in 
such a way that the Microsurface class can operate with both through a common interface. The random 
walks required to implement unbiased multiple scattering from the microsurface are implemented only 
once (each for eval() and sample()), regardless of NDF and microfacet material type. 

One benefit of this approach is that the Microsurface class is an operator that takes as input an NDF and 
a microsurface BSDF and outputs (implements) a new BSDF. This avoids having to specialize the random 
walks for each material type. For example, instead of implementing a RoughDielectricBSDF material, one 
simply passes the desired NDF and microsurface BSDF to the Microsuface class: 


DielectricBSDF micro_brdf; 
StudentTNDF ndf(&micro_brdf, rough_x, rough_y, gamma); 
Microsurface brdf(&ndf ); 


and then calls eval() and sample() with the two directions and iors on either side of the medium. Because 
the Microsurface is itself a BSDF, it can be passed into another Microsurface to model biscale roughness, 
etc. We have briefly explored this behaviour in the supplemental images, discussed in the next section. 

Our code is a direct extension of the prior implementation of Heitz et al. (2016) for heightfields, and 
the major differences are in the combination of support for singular and general microfacet BSDFs in the 
same random walk together with the joint sampling of new heights and microfacet orientations (which is 
required by the Null algorithm). 
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(a) Beckmann 0.4 ) Beckmann 1.0 (c) Beckmann 4.0 


(d) vMF 0.4 (e) vMF 1.0 (f) vMF 4.0 


Figure 7: Our Null-scattering formulation of rough surface scattering permits implementation of new 
microfacet BRDFs with full-sphere NDFs (in this case vMF) by implementing a single function: the NDF 
distribution itself. Both models are similar for low roughness, but we see the porous model transition 
naturally into high roughness regimes without the unnatural sheen appearance that Beckmann suffers as 
all of the microfacets are constrained to the upper hemisphere. 


5 Results 


5.1 Biscale Roughness 


In the supplemental images directory we provide a series of renders comparing biscale rough gold with 
Beckmann NDFs and identical roughnesses at each scale to a single-scale Beckmann roughness with an 
appropriately-scaled roughness to attempt an approximate match. We note that the biscale results have a 
less fuzzy/sheeny appearance, and that the difference is visible for surprisingly small roughness values. 


5.2 Student-T Conductors and Diffuse Surfaces with Multiple Scattering 


With an exact vNDF sampling for the Student-I NDF we computed (possibly the first) renders of rough 
conductors and diffuse surfaces with full multiple scattering. See the supplemental images folder to browse 
through the results over roughness and shape parameter and to compare the additional energy gained by 
including multiple scattering. 


5.3 Porous Rough Gold using the VMF NDF 


We implemented a NullINDF implementation of a spherical Gaussian (von-Mises-Fischer/vMF) NDF with 
roughness parameter a chosen to match Beckmann at low values, by implementing D as: 


virtual double D(const Vector3& wm) const 
{ 
const double u = wm.z; 
return exp((2.0 » (-1.0 + u)) / (m-roughness»m-roughness )); 


ys 


No additional code is required and we can now render porous rough gold materials, etc. by passing 
whatever microfacet BSDF we like into the NDF class and then into the Microsurface class. We compare 
rough gold renders for three roughness values {0.4, 1.0, 4.0} for Beckmann and vMF NDFs in Figure 7. 
Note how the Beckmann NDF unnaturally forces most facets to the vertical orientation for very high 
roughness values creating a bright sheen look, whereas the full-sphere NDF appears rougher and darker 
in a natural way as light bounces around in porous regions many times, being repeatedly absorbed. 
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